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ABSTRACT 

Some accreting neutron stars and young stars show unexplained episodic flares in the form 
of quasi-periodic oscillations or recurrent outbursts. In a series of two papers we present new 
work on an instability that can lead to episodic outbursts when the accretion disc is trun- 
cated by the star's strong magnetic field close to the corotation radius (where the Keplerian 
frequency matches the star's rotational frequency). In this paper we outline the physics of 
the instability and use a simple parameterization of the disc-field interaction to explore the 
instability numerically, which we show can lead to repeated bursts of accretion as well as 
steady-state solutions, as first suggested by Sunyaev and Shakura. The cycle time of these 
bursts increases with decreasing accretion rate. These solutions show that the usually assumed 
'propeller' state, in which mass is ejected from the system, need not occur even at very low 
accretion rates. 
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1 INTRODUCTION 

The interaction between a strong stellar magnetic field and an ac- 
cretion disc can affect both the evolution and observational proper- 
ties of the star. Close to the star the field is strong enough that the 
accretion disc is truncated, and mass is channelled along field lines 
to accrete on to the star's surface. At the inner edge of the truncated 
disc, the field and disc interact directly over some finite region, al- 
lowing angular momentum exchange from the differential rotation 
between the Keplerian accretion disc and the star. 

Angular momentum exchange between the field and the disc 
leads to two different states that can exist for a disc truncated by a 
magnetic field. The distinction depends on the position of the trun- 
cation radius relative to the corotation radius, r c = (GM./fi*) 1 ' 8 
(where M* and fi» are respectively the mass and spin frequency 
of the star), the radius at which the Keplerian frequency in the disc 
equals the star's rotational frequency. If the disc is truncated in- 
side r c then the field-disc interaction extracts angular momentum 
from the disc and accretion can proceed. If on the other hand the 
disc is truncated outside r c , the star-field interaction will create a 
centrifugal barrier that inhibits accretion. This is usually called the 
'propeller regime', under the ass umption that most of the ma ss in 
the disc is expelled as an outflow falarionov & Sunvaevll 19751) . 

Accreting stars with strong magnetic fields such as T Tauri 
stars, and X-ray millisecond pulsars show a large degree of vari- 
ability in luminosity (corresponding to changes in accretion rate), 
which may be ascribable to magnetic activity. For example, the pro- 
tostar EX Lupi (the prototype of the 'EXor' class), a TTauri star, 
increases a nd decreases in brightness by several magnitudes every 
1-5 years (Herbig|[2003). At much higher energies, a 1 Hz quasi- 
periodic oscillation (QPO) in accreting millisecond pulsar SAX 



11808.8-3658 has been observed during the decay phase of sev- 
eral outbursts dPatruno et ai]|2009h . The time-scale and magnitude 
of the variability in both sources suggest changes in accretion rate 
in the inner regions of the accretion disc, where it interacts with the 
star's magnetic field. 

In this paper we rev isit a disc instab ility first suggested in 
ISunvaev & Shakural d 19771) and developed in lSpruit & Taan] d 19931) 
(hereafter ST93), which can lead to episodic bursts of accretion. 
The instability arises when the magnetic field truncates the disc 
near the corotation radius. The magnetic field initially truncates the 
disc outside but close to the corotation radius, thus transferring an- 
gular momentum from the star to the disc and inhibiting gas from 
accreting on to the star (the propeller state). However, close to r c , 
the energy and angular momentum transferred by the field to the 
gas will not be enough to unbind much of the disc mass from the 
system and drive an outflow. Inste ad, the interaction with the mag- 
netic field will prevent accretion dSunvaev & Shakurj |l977l) . As 
gas in the inner regions of the disc piles up, the local gas pressure 
increases, forcing the inner edge of the disc to move inwards until 
it crosses r c . When the inner region of the disc cross inside r c , the 
centrifugal barrier preventing accretion disappears (since now the 
differential rotation between star and disc has changed sign) and 
the accumulated reservoir of gas is accreted on to the star. Once the 
reservoir has been accreted, the accretion rate through the disc's 
inner edge decreases, and the disc will again move outside r c , al- 
lowing another cycle to start. 

We study this process by following the time evolution of a 
thin axisymmetric viscous disc, with a paramaterization of the in- 
teraction between the disc and the magnetic field both inside and 
outside r c . This approach allows us to investigate the behaviour of 
the disc on time-scales much longer than the rotation period of the 



2 C. R. D'Angelo & H. C. Spruit 



star. Long time-scales are important since the instability evolves on 
viscous rather than dynamical time-scales of the disc. We are able 
to reduce the uncertainties in the detailed MHD interaction between 
the field and the disc to two free (but constrained) parameters. Us- 
ing this description we can then investigate the physical conditions 
for which the instability develops. 

In this paper we describe in detail the physics that can lead 
to episodic bursts of accretion and give a brief overview of the ob- 
served oscillations. In a later paper we will explore the range of 
outbursts seen in our simulations in more detail, and discuss their 
prospects for observability in specific stellar systems. 



2 MAGNETOSPHERE-DISC INTERACTIONS 

2.1 Interaction region between a disc and magnetic field 

We consider a star with a strong dipolar magnetic field surrounded 
by a thin Keplerian accretion disc. We assume that the dipole is 
aligned with both the star's spin axis and the spin axis of the disc, 
so that the system is axisymmetric. Near the surface of the star the 
magnetic field will truncate the disc, forcing gas into corotation 
with the star. This inner region (in which the gas dynamics is regu- 
lated by the magnetic field) is called the magnetosphere, and we de- 
fine the magnetospheric radius, r m as the radius at which the mag- 
netic fi eld is no longer stron g enough to force the disc into coro- 
tation jSprait&Taarr]! 19931) . Outside r m the magnetic field will 
penetrate the disc and become strongly coupled over some radial 
extent, which we call the interaction region, Ar. Beyond the inter- 
action region the disc and magnetic field are decoupled, so that the 
outer parts of the disc are not directly affected by the stellar mag- 
netic field. Figure [T] shows a schematic picture for the magnetic 
field configuration, with a closed magnetosphere close to the star, 
and a large region of opened field lines further out. 

In the interaction region, the differential rotation between the 
Keplerian disc and star shears the magnetic field, generating an az- 
imuthal component from the initially poloidal field. This in 
turn creates a magnetic stress which exerts a torque on the disc, 
transferring angular momentum between the disc and star. The 
torque per unit area exerted by the field on the disc is given by 
dtau/Ar — rS^z, where 

S z ^±^± (1) 

is the magnetic stress generated by the twisted field lines. The sign 
of the torque will depend on the location of the coupled disc re- 
gion relative to the corotation radius, r c = (GM* . If the 
coupling takes place inside r c the torque will extract angular mo- 
mentum from the disc, spinning down the disc (and spinning up 
the star), while if the coupling is outside r c the torque adds angular 
momentum to the disc, spinning it up (and spinning down the star). 

The radial extent of the interaction region has been a point 
of long-standing controversy in the study of accretion discs. 
In an early series of influential papers, Ghosh et al. (1977; 
lGhosh&Lamblll979al lbl) argued that the coupled region is large 
(Ar/r ^> 1), so that the magnetic field exerts a torque over a 
considerable fraction of the disc with a resulting large influence 
on the spin evolution of the star. However, the original mod el pro- 
posed by Ghosh & Lamb was shown to be inconsistent by IWand 
(1987). since the magnetic pressure they derived from field wind- 
ing far from r c is high enough to completely disrupt the majority 
of the disc. 

More recent analytical and numerical work has shown that the 



interaction region is likely m uch smaller, and much of the disc is 
disconnected from the star (seeJJzdensky 2004] for a recent review). 
This comes about from the fact that in force-free regions (where the 
magnetic pressure dominates over the gas pressure) as are likely to 
exist above an accreti on disc, field lines will tend to open u p as 
the twisting increases JaTvII 19851 ; iLvnden-Bell & Boilvll 19941) . As 
the disc and star rotate differentially, the increasing twist A<f> in 
the field line will only increase the B^ component to some maxi- 
mum B<f, ~ B z before the increased magnetic pressure above the 
disc causes the field lines to become inflated and eventually open, 
severing the connection between the di sc and star. Analytic studies 
of a sheared force-free magn etic field dAlvll 19851 ; Ivan Ballegooiienl 
1 1994 1 Uzdenskv et al. 2002) have shown that the B^ component 
will grow to a maximum twist angle A</> ~ tt before opening. 
The twist angle grows on the time-scale of the beat frequency 
= |f2* — Qk \~ , which is very short compared to the viscous time- 
scale in the disc except in a very small region around corotation. 

To prevent field lines from opening, they must be able to slip 
through the disc faster than the rate at which the field is being 
wound up. The rate at which the field can move through the disc 
is set by the effective diffusivity, r\, of the disc. Like the effective 
viscosity, v, that drives the transport of angular momentum, the 
effective diffusivity is also assumed to be driven by turbulent pro- 
cesses in the disc. Recent numerical studies of MRI (Magnetorota- 
tional Instability) turbulence (believed to be responsible for angular 
momentum transport in at least the inner regions of accretion discs) 
have tried to measure r\ directly. In these simulations, an external 
magnetic field is imposed on a shearing box simulation, and the ef- 
fective magnetic diffusivity is estimated as the flow becomes unsta- 
ble. The results suggest that the effective diffusivity and viscosity 
are of similar size, that is, t he effective magnetic Pr andtl number, 
Pr = ^/r/ is of order unity l lFromang~& Stone 2009). Such a large 
magnetic Prandtl number implies that for realistic disc parameters 
the magnetic field will not be able to slip through the majority of 
the disc fast enough enough to prevent fiel d lines from opening 
dLovelaceetaT1ll995l ; lUzdenskv et ai1l2002l) . Outside this region 
there will still be some coupling between the disc and the star as 
the gas moves from Keplerian to corotating orbits, but this estimate 
suggests that the actual extent of coupling is small (Ar/r < 1) 
regardless of where the disc is truncated relative to the corotation 
radius. 

Once the field lines are opened, there may be some re- 
connection across the region above the disc between open mag- 
netic field lines (e.g. lAlv & Kuijpersl Il990l ; iGoodson et all 1 19971 ; 
lUzdenskv et~afl l2002). The effective size of the interaction region 
would then depend on the efficiency of reconnection, and could also 
then become time-dependent (although likely on time-scales of or- 
der the dynamical time, which is much shorter than the viscous evo- 
lution time-scale). The opening and reconnection of field lines has 
also been suggested as a pos sible launching mecha n ism for strong 
disc winds and a jet (e.g . lAlv & Kuiipers] Il990l ; lHavashi et al.l 
ll996l ; lGoodson et aljfl997l) . This picture of a small interaction re- 
gion w ith some reconnection was first proposed bv lLovelace et al.l 
(1995), and has been supported by 2 and 3D si mulations of ac- 
cretion discs interacting with a magnetic field (e.g.lMiller & Stone 
1 1997 1; lGoodsonetaT1ll997l ; lHavashi et ai] 1 19961 ; iRomanova etal 

hocm 

In summary, although the extent of the interaction region is un- 
certain (subject to uncertainties in the effective diffusivity of mag- 
netic field in the disc and its possible reconnection in the magneto- 
sphere, as well as the detailed interaction between the disc and field 
near the magnetosphere), numerical and analytic work suggests that 
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Figure 1. Global magnetic field configuration for a strongly magnetic star 
surrounded by an accretion disc. In this picture, the majority of the field 
exists in an open configuration, and the connected region between the field 
and the disc is very small. Adapted from Lovelace et al. (1995). 



it is sma ll. Except for very special geometries for the m agnetic field 
(such as Agapitou & Papaloizou 200fl lShu et al.1 19941) . the low ef- 
fective magnetic diffusivity in the disc will force the magnetic field 
into a largely open configuration, and the majority of the accretion 
disc will be decoupled from the star, in strong contrast to the pre- 
diction of the lGhosh & Lambl dl979al) model. 

The extent of the interaction region as well as the average 
magnitude of the B<$, component generated by the disc-field inter- 
action will depend on the detailed interaction between the disc and 
the field as the gas moves from Keplerian orbits to corotation with 
the star, as well as the frequency and magnitude of possible recon- 
nection events. In the present work we therefore assume that the 
time-averaged B^ component generated by field-line twisting will 
be some constant fraction of B z , so that B^/B z = rj < 1. We also 
assume that Ar/r is small (< 1) but leave it as a free parameter. 



2.2 Accretion and angular momentum transport 

In this paper we describe the evolution of an accretion disc in 
which the conditions at the inner boundary are changing in time. 
Before doing this, however, we review how the conditions at the 
inner boundary affect the angular momentum transport and density 
structure of a thin accretion disc. In the thin-disc limit the evolution 
equation for the surface density E can be written: 

<9E 3 d 1/2 d , 1/2 .,. 
m=rd-r [r dr (r ^ (2) 

where v is the effective viscosity in the disc that enables angular 
momentum transport. In a steady state (in which the accretion rate 
is constant throughout the disc), the general solution for i^E is given 
by: 



i/E = 



7r ( ^ \ r / 



1/2 



(3) 



where is the inner edge of the disc, mis the accretion rate and /? is 
a dimensionless measure of t he angular momentum flux t hrough the 
disc p er unit mass accreted dPopham & Naravar]|l99ll : IPac zvnski 



All accretion discs have a boundary layer at their inner edge 
that connects the disc with either the surface of the star or the star's 
magnetosphere. In the boundary layer the gas must transition from 
Keplerian orbits to orbits corotating with the star in order to accrete. 
The structure of this boundary layer will determine the value of /3 
in {3). In the standard accretion scenario, that is, for accretion on to 
a slowly-rotating star or on to the star's magnetosphere inside the 
corotation radius, the gas in the boundary layer will be decelerated, 
meaning that there will be a maximum in the rotation profile, fi(r). 
At the maximum in Q(r), there is no longer an outward transfer 
of angular momentum from viscous torques, which in the thin-disc 
approximation w ill cause the surface density to decrease sharply, s o 
that /3 = 1 in OJ jPringle & Reeslll972l : [shakura & Sunvaevlll973l) . 
The maximum in Sl(r) effectively corresponds to the inner radius 
of the disc, since inside this radius gas is viscously decoupled from 
the rest of disc. The gas falling through the inner boundary of the 
disc will add its specific angular momentum (mrf n Q) to the star, 
spinning it up. 

However, there are in fact a wide range of solutions for the 
surface density profile of an accretion disc depending on the con- 
ditions imposed by the boundary layer, which in turn set the rate of 
angular momentum transport across the inner bou ndary of the disc . 
In a nonmagnetic star spin ning close to breakup dPaczvnskHll99ll : 
Poph am & Naravan|[l99lh . the angular momentum flux can be in- 
ward or outward, depending on the accretion history of the star. The 
dimensionless angular momentum flux /3 can in principle have any 
value less than 1 in this case. The top panel of Fig. [2] shows the 
steady-state surface density profile for a range of different values 
of /? from -1 to 1. 

ISunvaev & Shakural d 19771) studied a similar situation in 
which there is outward angular momentum transport in an accre- 
tion disc, and showed adding angular momentum at the inner edge 
of the accretion can in fact halt accretion altogether. The evolution 
of the disc in this case depends on the rate at which angular momen- 
tum is being injected at the inner edge of the disc compared to the 
rate at which it is carried outwards via viscous coupling. If angular 
momentum is injected into the inner boundary of the disc at exactly 
the same rate as viscous transport carries it outwards, then all ac- 
cretion on to the star will cease. For a steady state like this to exist, 
the outward angular momentum flux due to the magnetic torque at 
the inner edge of the disc has to be taken up at some larger distance. 
In a binary system, this sink of angular momentum can be the orbit 
of the companion star. If the disc is sufficiently large, the angu- 
lar momentum can also be taken up by the outer parts of the disc, 
while the inner parts of the disc are close to a steady state. The in- 
ner edge of the disc then slowly moves outward under the influence 
of the angular momentum flux. The surface density distribution in 
this case can be found from <[3j by taking the limit m — > 0, while 
letting f3 — !> — oo (noting that it measures the angular flux per unit 
accreted mass). This yields: 



r . 



(4) 



disc pi 
ll99l|). 



where /(n) is a measure of the torque exterted at the inner edge 
of the disc. The bottom panel of Fig. 2 shows the surface density, 
scaled to the value of /(n), for two instances of © with different 

valu es of rj n . 

ISunvaev & Shakural d 19771) refer to this solution as a 'dead 
disc', since there is no accretion on to the star. In this paper we 
call non-accreting discs without large outflows 'quiescent discs', to 
avoid confusion with 'dead zones' thought to be present in proto- 
stellar discs (regions in which there is insufficient ionization to 
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Figure 2. Surface density vT, of a thin disc as a function of distance from 
the corotation radius r c , for a steady, thin viscous disc. Top: steady accre- 
tion at a fixed accretion rate m, for inner edge of the disc at corotation. j3 
measures the angular momentum flux, (3 = 1 corresponding to the stan- 
dard case of accertion on to a slowly rotation object. For ft < the angular 
mometum flux is outward (spindown of the star). Bottom: 'quiescent disc' 
solutions with rh = and a steady outward angular momentum flux due 
to a torque / applied at the inner edge. The two curves show solutions for 
ri„/r c = 2 and 4. 



drive angular momentum transport via MRI but are too hot for ef- 
fici ent angular mo mentum transport via gravitational instabilities; 
e.g.lGammie 199(J). These quiescent discs play a role in the cyclic 
solutions discussed in Section[3] In these solutions accreting phases 
are separated by long intervals in which the inner disc is close to 
the quiescent state described by 10. 



2.3 Evolution of a disc truncated inside the corotation radius 

When the accretion disc is truncated by a magnetic field inside the 
corotation radius, the standard /3 — 1 case applies for a steady- 
state solution. The location of the inner edge of the disc ri n will be 
determined by the interaction between the disc and magnetic field, 
and change with changing conditions at the inner edge (such as the 
accretion rate on to the star). Here we estimate the location of ri n , 
and use it to show how the inner boundary of the disc will change 
in a non-steady disc. 

We define the inner edge of the disc as the point at which ma- 
terial in the disc is forced into corotation with the star. We use the 



azimuthal equation of motion for gas at the magnetospheric radius 
to obtain an estimate for ri n in a disc (see, e.g. ST93): 



2ir'S-^-(rv ( p) - ^--^-(rv^,) + 2irrS z<t> = 
at r or 



(5) 



where rhi n = — 2irrY,v r is the accretion rate through the inner 
edge of the disc. {5} neglects viscous angular momentum transport 
through the inner regions of the disc, under the assumption that it 
will be much smaller than angular momentum transport from the 
magnetic field. Using v^, — fi„r (since at n n the gas corotates 
with the star), and assuming a steady-state solution (d/dt = 0), 
(|5) becomes: 



rinS. 



4tt 



(6) 



where S z <t, is the magnetic stress from the coupling between the 
disc and star (introduced in Section [ZT). As long as the wind-up 
time for the field is shorter than the rate at which n n is changing, 
Btf, I B z will be roughly constant, so we make the assumption that 
Btj, = r/Bz, where r\ < 1 and is constant. 

For a dipole field aligned with the star's axis of rotation (B z = 
/i/r 3 , where /i = BsRt is the star's magnetic dipole moment), ^ 
can be re-written: 



( w 2 



77li r 



1/5 



(7) 



For 7] = 0.1, this estimate gives a value for n n about 40% smaller 
than the simple estimate found by equating the magnetic pres- 
sure from the field (B 2 /8ir) to the ram pressure fr om spherically- 
symm etric gas in free-fall on to the star (e.g. IPringle & Reesl 
|l97i). 

The derivation for n n above holds for steady accretion. For 
the problem studied here the position of the inner edge (set by the 
location of the magnetosphere) will change in time, which requires 
a minor reinterpretation of {7}. If n n is moving in time, the mass 
flux m co in the reference frame comoving with ri n differs from the 
mass flux, rh, measured in a fixed frame: 



= rh + 2nrT,r- ln 



(8) 



where fi n is the time derivative of n n . 

Since the torque between the magnetosphere and the disc acts 
at the inner edge, the mass flux entering the magnetosphere (used 
in ifTJ) is given by m co , not rh. As before, rh itself is given in terms 
of the surface density by the usual thin disc expression: 



3r 



1/2 



d_ 

dr 



.1/2 



(9) 



2.4 



Evolution of a disc truncated outside the corotation 
radius 



If the star is spinning fast enough, the magnetic field can truncate 
the disc outside r c . In this case the interaction with the magnetic 
field will add angular momentum to the disc, creating a centrifu- 
gal barri^rjhatijflubitsaccretion. This scenario was first described 
bv llllarionov & Sunvaevl i 1975b and is often termed the 'propeller' 
regime, under the assumption that the interaction with the magnetic 
field will expel th e disc at n n as an outflow v ia the 'magnetic sling- 
shot' mechanism ( Bland ford & Pavnell 19821) . 

However, in order for the gas to be ejected from the sys- 
tem, it must be accelerated to at least the escape speed (v e3C = 
\l GM,/2r). At the inner edge of the interaction region the gas is 
brought into corotation with the star, where v c — Sl,r. If this is 
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less than the escape speed, the majority of the gas will not be accel- 
erated enough to be expelled. Setting v esc = v c — \J GM*/r^r 
implies that for n n < 1.26r c most of the gas will not be expelled. 

Part of the disc could still be expelled in an outflow, but while 
the majority of the gas remains confined in the disc, the disc can 
act as an efficient sink for angular momentum from the star and ac- 
cretion can effectively be halted. The open field lines at larger radii 
could launch a disc wind which would provide an additional sink 
for angul ar momentum and so mewhat change the structure of the 
disc (e.g. lMatt & Pudritz 2005). Numerical studies of the field-disc 
interaction, for example, find that r econnection across fi eld lines 
can lead to intermittent accretion (e.g. Goo dson et al J 19971 see also 
Section[6](. However, models of disc winds typically include mass 
loss rate as a parameter of the problem, so that the amount of mass 
actually lost to the wind is uncertain. In this paper we make the as- 
sumption that the disc becomes quiescent, that is, for n n > r c no 
accretion or outflows occur. The steady-state disc solution is then 
given by I0. 

In the next section we will derive /(ri n ), the boundary con- 
dition for the surface density at the inner edge of a quiescent disc. 
Like for cases when n n < r c , we want to study non-steady-state 
solutions in which n n moves in time. As in the steady-state case, 
to derive n n we consider the difference in accretion rate at n n in a 
fixed frame and in a frame comoving with ri n . Since for a quiescent 
disc no matter is being accreted on to the star, m co = 0, so that {8} 
can be written: 



27rrEf in = Srll 2 ^ (r 1/2 z/E 



(10) 



Together with l[2}, a viscosity prescription and condition for the 
outer boundary, we can use the results from this section and the 
previous one to study the time-dependent behaviour of an accretion 
disc interacting with a magnetic field. 



3 CYCLIC ACCRETION 

The existence of quiescent disc solutions can naturally lead to 
bursts of accretion. Since there is very little accretion on to the star 
or outflow, if mass continues to accrete from larger radii it will 
pile up in the inner regions in the disc until the gas pressure is high 
enough to overcome the centrifugal barrier from the magnetic field- 
disc interaction and accretion can proceed. Once the reservoir has 
been emptied the inner edge of the disc will move back outside the 
corotation radius and the reservoir will start to build up again. 

In Sections |2.3I and [2~4l we showed how the inner radius of a 
thin viscous accretion disc will evolve inside and outside corota- 
tion. To study the time-dependent evolution of a disc, we must con- 
nect these two states as the inner edge of the disc passes through 
the corotation radius. We also require a description for /(nn), the 
inner boundary condition for the disc truncated outside r c . 



3.1 Surface density profile for n n > r c 

When the interaction region is outside r c , the star is rotating faster 
than the Keplerian disc and the magnetic field lines lead the disc, 
adding angular momentum to the material in the inner regions. As 
discussed in Section |2~T1 the torque per unit area exerted on the 
disc will be {S^, z )r, so that the torque exerted across the entire 
interaction region (assuming it is small) is approximately: 



where the extra factor 2 comes from coupling to both sides of the 
disc. 

As argued in the previous section, if the disc is truncated close 
to but outside r c , the majority of the gas in the interaction region 
will not be expelled in an outflow. Instead, the angular momen- 
tum from the magnetic field is transferred outwards to the rest of 
the disc. We can derive a relationship between the position of and 
surface density at the inner edge of non-interacting disc from the 
conservation of angular momentum across the interaction region. 

Since the interaction region is small we do not consider its 
density profile explicitly, focusing instead on its influence on the 
non-interacting disc. We therefore define n n as the point in the disc 
just outside the interaction region, where there is no magnetic cou- 
pling between the disc and the star. Across the interaction region 
the density in the disc decreases sharply (since the gas is forced 
into nearly corotating orbits with the star). We make the simplify- 
ing assumption that none of the mass in the disc escapes, either into 
an outflow or through the magnetosphere on to the star. The inner 
edge of the interaction region, r- m — Ar, is therefore defined as the 
point at which the surface density drops to zero. 

To determine E at n n we consider the angular momentum flux 
across Ar when n n > r c . The flux of angular momentum must 
be continuous across Ar, meaning that the viscous angular mo- 
mentum transport outside Ar must balance the angular momentum 
flux added by the magnetic field across the interaction region. This 
balance is written: 



mr 2 n - 2irr(vT,) + r 2 n' = 

rhr 2 Q - 27rr(^E)~r 2 fi' + 



(12) 



, + Ar 



47rr Sz^dr. 



In this equation, v and E* are the viscosity and surface den- 
sity inside ( _ ) and outside ( + ) Ar, m = 27rr(Ew r ) ± is the mass 
flux through Ar (where v r is the radial velocity of the gas) and fl 
is the orbital frequency at n n . The first term on either side of the 
equation denotes the advection of angular momentum across ri n , 
while the second is the angular momentum transported by viscous 
stresses. The final term on the right hand side is the angular mo- 
mentum added by the magnetic field to the coupled region of the 
disc. The first term on both sides cancel (to enforce conservation 
of mass across Ar), and we make the further assumption that in 
the interaction region most of the angular momentum is transported 
through external magnetic torques rather than viscous stress, so that 
(j/E)~ <C (i/E) + . For a small interaction region, the last term in 
d 1 2b can be re-written: 



4irrS z4 ,dr w 47rArr in (5 z ^). 



(13) 



J12t can then be re-written to yield the surface density at 

rin for r > r c : 

2(5^)Ar 



(14) 



r ~ 47r(S , z )r i nArz, 



(11) 



As predicted in Section [2~4l d 14b shows that the surface den- 
sity at n n will be large, a consequence of the torque being applied 
by the disc-magnetic field coupling iSunvaev & Shakuri] 1 19771 : 
iPopham & Naravarj 1 199 ll : |Paczvnski|[l99lh . JT4b corresponds to 
the function /(ri n ) introduced in Section [2~2l for n n > r c , that 
is, the boundary condition at the inner edge of the disc. In a time- 
dependent system, as gas accretes from larger radii (via viscous 
torques) it will pile up near n n and the increased gas pressure will 
push the inner edge of the disc further inwards towards r c . 
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3.2 Transition region 



3.3 Physical constraints on Ar and Ar 2 



When the inner edge r- m is well inside r c , conditions at the inner 
edge are the standard ones for accretion of a thin disc on a slowly 
rotating object: 



E(r in ) = 0, 



(15) 



while the time-dependent position of the inner edge is determined 
by Q: 



1/5 



(16) 



where m co is the mass flux in a frame comoving with r- m as dis- 
cussed above. 

When the inner edge is outside the corotation radius, the mag- 
netosphere does not accrete: 



m co = 0, 



(17) 



while the surface density at n n is determined by a magnetic torque, 
as discussed above. With the Keplerian value for fi(ri ft ) and as- 
suming a dipolar magnetic field, the results of Section [3~T1 can be 
re-written: 



Ar 



37T(GAf») 1 /2 r 9/2' 



(18) 



To connect these two limiting cases, we assume that the ef- 
fect of the interaction processes is equivalent to a smooth transition 
in the conditions. This is valid since the time-scales we are inter- 
ested in are much longer than the orbital time-scale on which the 
conditions of the transition region between disc and magnetosphere 
vary. The assumption is thus that the effect of the fast processes in 
the transition region can be represented by averages. The mass flux 
on to the magnetosphere is therefore taken to vary smoothly from 
for Tin well outside corotation to the value in dl6t valid well inside: 



(19) 



where m + is given by l !16t . For the connecting function y m we take 
a simple function that varies from to 1 across the transition: 



tanh 



Ar 2 



(20) 



where Ar 2 is the nominal width of the disc-magnetosphere transi- 
tion and a parameter of the problem. 

Similarly the surface density at the inner edge makes a smooth 
transition from its value in d 1 8b to 0: 



E in = V T,TP 



(21) 



In this paper we treat Ar and Ar 2 as free parameters. However, a 
lower limit on both parameters can be set by considering the stabil- 
ity of the inner regions of the disc to the interchange instability. In 
the quiescent disc, the low-density magnetosphere must support the 
high-density disc against infall. This configuration will be unsta- 
ble to interchange instability (the analog of the Kelvin-Helmholtz 
instability), unless the surface density gradient in the interaction 
region is shallow enough to suppress it. This sets a limit on the 
minimum width of the interaction region, Ar, where the density 
gradient falls from its maximum (at ri n ) to close to zero in the mag- 
netosphere. 

This instability also sets a limit on the minimum width of 
Ar 2 , the transition length over which the disc moves from a non- 
accreting quiescent disc to one in which there is accretion through 
the inner boundary. As ri n moves closer to r c the width of the in- 
teraction region preventing accretion (i.e. where the field lines are 
adding angular momentum to the disc) decreases. When the width 
of the interaction region outside r c becomes smaller than is stable 
against the interchange instability, accretion through the magneto- 
sphere will begin. Ar 2 must therefore be larger or equal to this 
value, that is, at this minimum distance from r c accretion onto the 
star will take plac e. 

ISpruit et alj (l995) studied the stability of a disc interacting 
with a magnetic field to interchange instabilities, and derived the 
following linear stability criterion: 



B r B z d 
2ttE dr 



In 



E 



> 2 



(23) 



Assuming that B r ~ in our formulation this inequality 
becomes: 



3a 



1 + tanh (3?) 



I > 2 I 1 | 



3/2 



(24) 



For a = 0.1 and assuming H/r is in the range 0.07-0.1, 
the range of Ar/r = [0.05,0.1] will satisfy this inequality for 
Ar 2 /r = [0.01, 0.02]. In this inequality larger values of Ar cor- 
respond to smaller possible values for Ar 2 , since larger Ar corre- 
spond to smaller maximum E(ri n ) and hence shallower gradients. 
This instabi lity has recently been studi ed using 3D numerical sim- 
ulations (Kulkarni & Romanoval 120081) . who find numeric a lly ap - 
proximately the same criterion for stability as lSpruit et alj dl995h - 
The shaded regions of Figs. [7] and [8] show the values for Ar 2 and 
Ar that are unstable to the instability studied in this paper. The sim- 
ple analysis of this section suggests that at least part of the shaded 
sections in Figs.[7]and[8]will be stable against the interchange insta- 
bility, so that the larger magnetosphere-disc instability could occur. 



where the connecting function j/e is: 



Vs = i [l + tanh ( — 



Ar 



(22) 



All the uncertainties in the transition region are thus subsumed 
in the parameters Ar and Ar 2 . In Section[5]we study the effect of 
these uncertainties with a parameter survey. The effective widths 
of the transition of magnetospheric accretion rate and inner-edge 
surface density need not be the same, and we in fact find that the 
difference between Ar and Ar 2 is important for the form of the 
resulting accretion cycles. 



4 NUMERICAL IMPLEMENTATION 

4.1 Disc equation and viscosity prescription 

To study the surface density evolution of an accretion disc inter- 
acting with a magnetic field as outlined in the previous section, 
we use a time-dependent numerical simulation of a diffusive ac- 
cretion disc. Our assumption that the interaction region is small 
(Ar/r < 1) means that rather than calculate the disc behaviour in 
the interaction region explicitly we can instead use the physics of 
the interaction region to derive boundary conditions for the inner 
edge of the non-interacting disc. 
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We assume that the accretion disc (outside the interaction re- 
gion) can be treated in the thin-disc limit, so that the evolution equa- 
tion for the surface density E is given by 0. We assume that the 
viscosity in the disc follows a power-law dependence, so that: 



for 



where v Q = a(GM,) 1/2 (H/R Y and 7 



standard a-viscosity prescription dShakura & 



(25) 

0.5 following the 
To 



U.o to liowing 
Sunvaev|[l973h . 



evolve ((2) in time, we require boundary conditions at n n and r out , 
plus an additional equation to describe the movement of the in- 
ner edge of the disc, fin- We set the outer boundary by defining 
the mass accretion rate through the outer edge of the disc (rh), 
which we vary as a parameter of the problem. This defines the time- 
averaged mass accretion rate in the disc. The surface density at the 
inner edge of the disc is given by J2U : 



E(r in ) = 



Ar 



(26) 



We calculate the displacement of the inner boundary using the 
results of Sections 12.31 and [2~4l by considering the difference be- 
tween the total mass flux at n n in a fixed and comoving frame of 
reference: 

m co = m + 27rrEri n , (27) 
where m co is given by H9i . This expression can be re-written: 



67rr j L / 2 ^(p'Eri n ) = -27rnnE(nn)nn + 



(28) 



1 — tanh 



Ar 2 



TjfJ, 



Taken together, (O, d25t . ( 126b . J28b and an outer boundary 
condition describe the time-dependent evolution of an accretion 
disc. 

4.2 Steady-State solution 

From the results of the previous sections, we can calculate the 
steady-state solutions for a given rh, the average mass accretion 
rate. For certain values of rh, Ar and Ar 2 , this equilibrium is un- 
stable, leading to oscillations in ri n and corresponding accretion 
bursts. 

In a steady-state, the accretion rate is constant throughout the 
disc, i.e. rh co — rh: 



1 — tanh 



r/fi 



(29) 



4firf n ' 

yields the inner radius of the disc in 



Implicitly solving j29\ for ft 
a steady-state solution. 

The general steady-state surface density profile was calculated 
in Section [2~2l and is given by 101 with an additional term since 
rh 7^ in the disc. The function /(nn) is given by equation J 1 8t >. 
The steady-state surface density profile will thus be: 



m 
3^ 



1/2 



(30) 



+ 



r/fj, Ar 



[1+ tanh 



6ivr? n (GMry/ 2 

The numerical simulations described in the following sections 
of the evolution of a viscous accretion disc show that the equi- 
librium solution given by j29\ and i30\ can become unstable to 
episodic bursts of accretion by the process outlined in Section[3] 



4.3 Numerical setup 

To follow the time-dependent evolution of a viscous accretion disc 
interacting with a magnetic field we use a ID numerical simulation, 
first making a series of mathematical transformations. 

The power-law prescription for the viscosity, l |25t , allows us 
to define a new function, u, for convenience: 



Er 



l/2+ 7 



(3D 



To make our results more readily applicable to different mag- 
netic stars (e.g. neutron stars, magnetic white dwarves and proto- 
stars), we adopt scale-free coordinates. The instability studied in 
this paper varies on viscous time-scales of the inner disc, which are 
in general much shorter than the time-scale over which the transfer 
of angular momentum between the star and the disc can substan- 
tially change the star's rotation period. A constant rotation period 
implies that a constant corotation radius, making it a natural choice 
for scaling our variables. We thus scale the radial coordinate to the 
corotation radius, and the time in terms of the viscous time-scale 
(r 2 /v) at the corotation radius. Further, since we are most inter- 
ested in the behaviour of inner regions of the disc, we adopt a co- 
ordinate system comoving with ri n : 



:t = t 



.2-7 ' 



(32) 



Dropping the prime notation, the surface density evolution equation 
in the new coordinate system then becomes: 



2a = 3r 7-v2 jL 

dt dr 



.1/2 



dv 
dr 



du 

+ r in — , 
or 



with the boundary condition at n n given by: 
■ 2 Ar 



u(ft n ) 



3-K(GM t y/ 2 vor* r in ln 



tanh 



+ 1 



The evolution of the inner edge of the disc is given by: 



tanh 



Ar 2 



TjfJ, 



-11/2+7 



167rO*i^ r r 3/2 ^Oin) 



(33) 



•(34) 



(35) 



3r^ du I 
""(n n ) dr \r ia ' 

Finally, to increase the resolution at the inner edge of the disc 
we make a further coordinate transformation to an exponentially 
scaled grid: 



In 



Tout 



1 



(36) 



where a is a scaling factor to set the clustering of grid points around 

n n . 

We calculate the second-order discretization of the spatial 
derivatives on an equally-spaced grid in x. To evolve the result- 
ing system of equations in time requires an algorithm suitable for 
stiff equations. This is necessary to follow the evolution of the in- 
ner boundary, l |28l l. When rjn S> r c , j2S\ reduces to a differential 
equation that is first order in time. However, for n n r C) E(n„) 
becomes very small, and the equation essentially becomes time- 
independent. We have formulated the problem so that E(n n ) stays 
small but non-zero for all values of nn (so that the solutions is con- 
tinuous at all values of nn), but its small value inside r c means that 
the differential equation is stiff (since the evolution equation for nn 
in l |28t contains terms of very different sizes). To perform the time 
evolution, we therefore use the semi-implicit extrapolation method 
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dPress etal.Ul992h .P. 724), which is second-order accurate in time 
and suitable for stiff equations. 

Since the grid comoves with the inner radius, the outer bound- 
ary of our disc also moves. We set the accretion rate at the outer 
boundary to be fixed in the moving coordinate system, so that it 
changes slightly as the outer boundary moves. The effect is neg- 
ligible as long as the disc is large enough that the outer parts of 
the disc are unaffected by the changing inner boundary condition, 
which we confirm by varying the position of the outer boundary of 
the disc. 

The solutions are sensitive to the changing conditions at the 
inner boundary of the disc. To confirm that our results are robust 
for the grid we have chosen, we varied the various numerical pa- 
rameters of the problem: grid resolution, the exponential stretch 
parameter a at the inner boundary (see l |36t ) and the fractional ac- 
curacy of the solution computed by the semi-implicit extrapolation 
method (which sets the maximum possible timestep). 



5 RESULTS 

Our primary goal in this paper is to study the conditions for which 
the disc is unstable to episodic outbursts. To do this we follow the 
evolution of an accretion disc in which the mean mass accretion 
rate, m is a parameter of the problem by setting rh as the accretion 
rate through the disc's outer boundary. The other system parame- 
ters of the problem are the stellar mass, M*, frequency, and 
magnetic moment, fj,, The interaction between the magnetic field 
and the disc introduces three additional parameters: r\ = B^/Bz, 
the fractional width of the interaction region Ar/r, and the length 
scale Ar2/r over which the inner edge of the disc moves from a 
non-accreting to accreting state. Finally, our description of the vis- 
cosity, J25K introduces three additional parameters: a, the aspect 
ratio of the disc, H/R (assumed constant), and 7, the radial power- 
law dependence of the viscosity. 

The problem has two scale invariances, which reduces the 
number of free parameters. As seen in {25}, a and H/R are de- 
generate. Additionally, the system parameters fi, M*, f2» and rh 
can be re- written as the ratio m/rh c , where rh c is the accretion 
rate in {7} that puts the magnetospheric radius at r c . This ratio is 
equivalent to the 'fastness parameter', fii n /fi* (where £l[ n is the 
Keplerian frequency at n n ) which is sometimes used to describe 
disc-magnetosphere interactions. 

For reference, our dimensionless parameter rh / m c can be ex- 
pressed in terms of physical parameters appropriate for protostellar 
systems: 



m 

m c 



2.3 x 10- 7 M©yr- 



f B s 



R, 



0.6M Q 



5/3 



7/3 



(37) 



V 200007 \2.ir q j V lda y, 

We assume that the time-averaged B^ component will be 
constant with radius in the coupled region, and set the parameter 
7/ = 0.1. For the viscosity, v = a(GKh) 1/2 (H / 'R) 2 r~< , we take 
a = 0.1 and H/R = 0.1 to calculate the magnitude of uq, and 
assume 7 = 0.5 everywhere in the disc. Varying a, H/R and 7 
will change the time-scale over which outbursts occur, but will not 
change the general character of our outburst solutions. 

This leaves three scale-free parameters in the problem: rh/rh c , 
Ar/r and Ar2/r. We vary each of these parameters to explore 
the range of unstable solutions. For small values of Ar/r (~ 0.1) 




0.05 



Time ft 



Figure 3. Growth of instability from steady-state solution, (29) and (30) , 

for m/riic = 1, Ar/r = 0.05, and Ari/r = 0.014. The inner radius 
(solid curve) evolves around its steady-state value (dashed horizontal line), 
causing the net accretion rate on to the star to change as well (dashed curve). 



and Ar-i/r (~ 0.01), and rh/rh c < 1, the position of the inner 
boundary quickly becomes unstable and begins oscillating. Since 
the position of r ln determines the mass accretion rate on to the star, 
(19) , the change in r- m leads to an accretion outburst. We use the 
steady-state solution (given by j29\ and OOt ) as an initial condition 
for all our simulations. 

Fig. [5] shows the growth of the instability for rh/rh c = 1, 
Ar/r — 0.05 and Ar%/r = 0.014. The solid curve shows the evo- 
lution in ri n , scaled to the corotation radius. The horizontal dashed 
line shows the steady-state value for r- m . The right-hand axis plots 
the accretion rate on to the star as a function of time (the dashed 
curve). The accretion rate is scaled to units of the steady-state ac- 
cretion rate, rh. The instability quickly grows out of the equilibrium 
solution, and saturates into steady oscillations. 

We observe a wide range of oscillatory solutions that span 
three orders of magnitude in frequency, depending on the values 
of rh/rh c , Ar/r and Ar-^/r. The shape of the accretion burst itself 
also changes dramatically depending on the system parameters. At 
large rh/rh c the bursts are quasi-sinusoidal oscillations, as in Fig. 
[3]and the bottom panel of Fig. [4] As the mean accretion rate is de- 
creased, the bursts take the shape of a relaxation oscillator, where 
the bursts are characterized by an initial sharp spike of accretion 
which then relaxes to a quasi-steady accretion rate for the duration 
of the burst, before abruptly turning off as the reservoir is emp- 
tied and ri n quickly moves well outside r c . During the outburst 
phase, higher frequency sub-oscillations are also sometimes seen 
with varying intensity. 

Figs. [4] and [5] show the evolution of ri n and accretion rate 
as we vary rh/rh c but the other parameters stay fixed. From bot- 
tom to top, the panels of Fig. [4] show the instability for rh / rh c 
= [0.095, 0.052, 0.031] (rh = [2.2,1.2,0.73] x lO^Meyr" 1 
for the parameters in {37)). At the highest mean accretion rate, 
ri n (the solid curve) oscillates with a high frequency around its 
steady-state value (dashed line), with corresponding bursts of ac- 
cretion on to the star (dashed curve). As rh/rh c is decreased, 
the accretion profile changes to much lower frequency outbursts, 
with long periods of quiescence as r- ln moves away from r c and 
accretion ceases completely. The high-frequency oscillation that 
dominates for rh/rh c = 0.095 is superimposed over the low- 
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Figure 4. Outburst profiles of r- ln and m for moderate values of rh/m c . From bottom to top, m/m c = [0.095,0.052, 0.031]. For adopted protostellar parameters 
this corresponds torn = [2.2, 1.2, 0.73] X 10 — 8 M0yr — x . The lines are the same as in Fig. [5] 



frequency accretion bursts for lower m/m c . Fig. [5] shows the 
continuation of Fig.|4]for m/m c = [0.019,0.0084,0.003,0.0022] 
(m = [4.5,1.9,0.95,0.38] x 10" 9 Moyr -1 ). The characteris- 
tic accretion burst profile essentially stays the same as m/m c is 
decreased, with sharp spikes at the beginning and end of an ac- 
cretion outburst. The overall amplitude of the outburst decreases 
only slightly with decreasing mean accretion rate. The initial spike 
decreases by about 20% as the mean accretion rate drops from 
rh/m c = 0.052 to rh/m c = 0.0022. The more significant effect is 
that the length of time between outbursts increases with decreas- 
ing m/rhc, since at low average accretion rates it takes longer to 
build enough mass to drive another outburst. The overall shape 
of the outburst is relatively insensitive to changing rh/m c , be- 
coming shorter as m/m c decreases. At the lowest accretion rate 
(3.8 x lO -lo M0yr -1 ; the top panel of Fig. [5}, the burst consists 
of only one sharp spike. As we have formulated the problem, the 
instability will persist down to arbitrarily low accretion rates. 

Changing the other parameters, Ar/r and Ar^/r, has a much 
stronger effect on the shape of the outburst than changing the mean 
accretion rate. Fig.[6]shows the outburst profiles for different values 
for Ar/r, setting rh/rh c = 0.04 and Ar2/r — 0.014. From the bot- 
tom to top, Ar/r = [0.03,0.05,0.07,0.09], which spans the unstable 
region of Ar/r for the adopted m/m c . For small Ar/r the insta- 
bility manifests itself as repeating short bursts of accretion, with 



comparatively long quiescent phases. As Ar/r increases, the fre- 
quency of the outburst decreases, and the duty cycle increases dra- 
matically. For very large Ar/r the outburst lasts about 200 times 
as long as for the minimum Ar/r but at lower accretion rate after 
the initial spike. The burst profile of the instability is thus sensitive 
to small changes in Ar/r, but the range in Ar/r over which the 
instability exists is quite small.We find a similar range of outburst 
profiles by changing Ar-2/r and keeping Ar/r fixed, except with 
the opposite trend: for large Ar2/V the instability manifests as a 
series of short spiky bursts, becoming longer as Ar2/V decreases. 

We next considered the parameter space in rh/rh c , Ar/r and 
Ar-^/r over which the instability occurs. We have briefly explored 
the effect of varying both Ar/r and Ar-^/r over a small range 
in m/m c and found that, although the outburst profile changes 
somewhat, the range over which Ar/r and Ar^/r produce un- 
stable solutions are independent. We therefore assume that Ar/r 
and Ar2/r vary independently of each other for all m/m c , and 
consider the range of the instability over the [m/m c , Ar/r] and 
[m/m c , Ar2/V] spaces separately. 

Fig.|7]shows the range of unstable solutions (shown as shaded 
regions) changing rh/rh c and Ar/r, but keeping Ar^/r fixed at 

0. 014. Although there is a small unstable branch around rh/rh c = 

1, in general as Ar/r increases, a lower rh/rh c is required before 
the instability sets in. 
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Figure 7. Parameter map of instability as a function of m/ rh c and width of 
interaction region Ar/r, with constant Ar2 = 0.014. The shaded regions 
denote unstable parameters. 
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Figure 8. Parameter map of instability as a function of rh/rh c and accretion 
transition length Ar2/r, with constant Ar2 = 0.014. The shaded regions 
denote unstable parameters. 



Fig. [8] shows the unstable solutions changing m/m c and 
Ar2/r but keeping Ar/r fixed at 0.05. The opposite trend from 
Fig.|7]is seen, with a larger range of unstable accretion rates. There 
is again a range of unstable solutions around m/m c = 1, although 



in this case the unstable region extends over the entire Ar-^jr pa- 
rameter space. The instability likely extends to smaller Ar2 /r, but 
we do not explore the region smaller than Art = 0.005 on physical 
grounds, since such a small transition length will likely be unsta- 
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ble to other instabilities like the interchange instability (see Section 
13.31 1. As with changing Ar/r, the outburst profile changes substan- 
tially over the small range of Ar2 /r in which the instability occurs. 



6 DISCUSSION 

In this paper we studied a disc instability first explored by 
ISunvaev & Shakural Jl977h and ST93, with a more physically mo- 
tivated and general formulation of the problem than was used in 
ST93. In particular, we have improved the description of the disc- 
field interaction when the disc is truncated outside corotation by 
deriving conditions for a 'quiescent' state, in which the angular 
momentum transferred from the star into the disc halts accretion 
altogether. In agreement with ST93, we observe a wide range of 
oscillatory behaviour, and the frequency range of individual out- 
bursts spans three orders of magnitude. 

The period of the cycle seen in Figs. [3}{5] varies from 0.02 
to 20t c , where tc is the nominal viscous time-scale at the coro- 
tation radius t c = r^/u(r c ). Though cycle times scale with i c , 
this is evidently not the only factor. As discussed in ST93, the vis- 
cous time-scale relevant for the cycle period depends on the size of 
the disc region involved. This depends itself on the cycle period, 
hence the period must be determined by additional factors. One of 
these is the mean accretion rate, but the physical conditions in the 



magnetosphere-disc interaction region have an equally important 
effect. 

From Figs.|7]and[8]it appears that there are two different kinds 
of instability. One of these operates in a narrow range of accretion 
rates, around the value where steady accretion would put the inner 
edge at corotation. The instability in this case is of the type shown 
in Fig. [3] an approximately sinusoidal modulation, characteristic 
for a weak form of instability. The inner edge of the disc oscillates 
about a mean value, but stays inside the width of the transition re- 
gion. The longer cycles in the upper parts of Figs. [7] and [8] are a 
strongly non-linear, relaxation type of oscillation. The inner edge 
is somewhat outside the transition region for much of the cycle with 
no accretion taking place (the 'quiescent' phase), and dips in for a 
brief episode of accretion befor e moving back out again. This is 
the kind of cycle envisaged by Sunyaev & Shakura ( 1977). Dur- 
ing the quiescent phase, the disc dSunvaev & Shakura l ll977h call 
it a 'dead disc') extracts angular momentum from the star by the 
magnetic interaction at its inner edge. These two forms of instabil- 
ity are merged into a continuum in ST93, as a result of the different 
(and less realistic) assumptions made there about the interaction be- 
tween disc and magnetosphere outside corotation. This difference 
also affects the dependence on the mean accretion rate. Whereas 
in ST93 cyclic behavior was found only in a limited range of ac- 
cretion rates, our results show that cycles can occur in principle at 
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arbitrarily low accretion rates, with steadily increasing cycle period 
and decreasing duty cycle of the accretion phase. 

Figs. [4] and [5] show that the radius of the inner edge of the 
disc does not move by more than 10% around corotation, even at 
the lowest mean accretion rates. For example in the case m/rh c = 
9.5 x 10 -2 of Fig. [4] the standard 'ram pressure' estimate would 
yield a much larger magnetosphere radius, about r m = 3.6 r c . The 
difference arises because in our cyclic accretion states the condi- 
tions in the inner disc are very different from those assumed in 
conventional estimates of r m ; the density in the inner disc, for ex- 
ample, is much higher. 

At ri n ^ l.lr c , the velocity difference between the magne- 
tosphere and the disc is only 5%, much less than the 40% which 
mass would need in order to escape from the system. 'Propeller- 
ing' of mass out of the system is thus unlikely to be effective. This 
does not exclude that some mass loss (powered by a magnetic wind 
from the disc or the interaction region around the inner edge of 
the disc) may also take place, but our results show that this is not 
a necessary consequence for a disc in what is traditionally called 
'propeller' regime. 

At sufficiently low accretions rates one would expect, how- 
ever, that propellering would also be a possible outcome: if the ro- 
tation rate of the star is high enough, matter could be ejected before 
it has the time to form a dense disc. The existence of a cyclic form 
of accretion at low accretion rates thus suggests that two different 
accretion states are possible, and that there would be a second pa- 
rameter determining which of the two is realised. This might simply 
be the history of the system. 

If a disc is initially absent and accretion is started, the den- 
sity will initially be low enough that ejection by propellering can 
prevent accretion altogether. The cataclysmic variable AE Aqr (e.g. 
IWvnn et alj|l997h is likely to be such a case. On the other hand, if 
a disc is initially in a high accretion state such that the inner edge is 
inside corotation, a subsequent decline to low accretion rates could 
lead to the cyclic accretion described here. Such a situation could 
be at work in the TTauri star EX Lupi (where the initial high ac- 
cretion phase has ended). It could also be appropriate for the X-ray 
millisecond pulsar, SAX J1808. 8-3658, whi ch has shown a 1-H z 
QPO in the decline phase of several outbursts IPatruno et alj (2009). 
The pile-up of mass at the magnetosphere will maintain the disc 
this state, and prevent propellering even when the mean accretion 
rate drops to very low values. 

The instability studied in this work has not yet been observed 
in numerical simulations, partly because most numerical simula- 
tions do not run for long enough to observe it, but mainly because 
most simulations have focused on either accreting or strong pro- 
peller cases. However, in virtually all numerical simulations out- 
flows and variability in the disc are observed, with an intensity that 
varies between different simulations. Gas pile-up at the inner edge 
of the disc is also observed, with the amount of pile-up tied to 
the effective diffusivity of magnetic field at the inner edge of the 
disc (e.g. ?). The process of closing and opening field lines pro- 
vides a source of mass to launch both a weakl y-collimated outflow 
(the disc wind) and a well-collimated jet (e .g lHavashi et alj[l996l ; 
Good son" et alJl997l ; lRomanova et alj2009f) . The whole cycle takes 
place on time-scales that can vary between the dynamical and vis- 
cous time-scales at the inner edge of the disc, but are generally of 
higher frequency than the disc instability studied in this paper. The 
inner edge of the disc also oscillates significantly (although it re- 
mains on average outside co rotation), from betwe en a few stellar 
radii jRomanova et alJl2009T) up to 30 stellar radii (Goods on et al.1 
1 1997b . Even if such variability is present, the instability studied 



in this paper can still occur provided the outflows/accretion bursts 
generated by field lines opening are not strong enough to fully 
empty the reservoir of matter accumulating just outside r c . 



7 CONCLUSIONS 

We have studied the accretion of a thin viscous disc on to a mag- 
netosphere of a magnetic star, under the influence of the magnetic 
torque it exerts on the disc. We focused in particular on cases with 
low accretion rates. For high accretion rates such that the inner 
edge ri n of the disc is inside the corotation radius, standard steady 
thin viscous disc solutions are recovered. However, when the inner 
edge is near corotation we find that the accretion becomes time- 
dependent, and takes the form of cycles consisting of alternating 
accreting and non-accreting ('quiescent') states. The period of this 
cycle varies from a small fraction of the characteristic viscous time 
scale in the inner disc, rf n /v, to a large multiple of it, depending 
on the mean accretion rate as well as on the precise conditions as- 
sumed at the magnetosphere. 

These cyclic accretion solutions continue to exist indefinitely 
with decreasing accretion rate. The cycle period increases, while 
the duty cycle of the accreting phase decreases with decreasing ac- 
cretion rate. In the quiescent phase after a burst of accretion, the 
inner edge of the disc moves outward, and mass starts piling up in 
the inner regions of the disc. In response, the inner edge eventually 
starts moving back in again and accretion picks up as n n crosses the 
corotation radius. This empties the inner regions of the disc, caus- 
ing the inner edge to move outward again. The cycle thus has the 
properties of a relaxation oscillator, as found before in ST93. The 
reservoir involved is the ma ss in the inner region of the disc. These 
results (as well as those of ISunva ev & Sha kural i 19771) and ST93) 
show that accretion without mass ejection can occur at accretion 
rates well inside what is usually called the 'propeller' regime. In- 
stead of the mass being ejected, the accreting mass can stay piled 
up at high surface density in the inner disc, just outside corotation. 
We have suggested that systems with very low accretion rates can 
be in either of these states. Propellering would occur when a disc 
is initially absent and mass transfer is first initiated (the case of AE 
Aqr for example), while a system with an accretion rate that drops 
from an initially high value would end in the cyclic accretion state 
described in this paper. This would apply to most cataclysmic vari- 
ables and X-ray binaries, as well as some TTauri stars. 
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